Apparatus and method for smoothing random event data obtained from a positron emission tomography scanner

ABSTRACT

A method and apparatus for smoothing random event data obtained from a Positron Emission Tomography (PET) scanner. The method includes obtaining initial random event data u(s, φ, t=0)=u 0 (s, φ), corresponding to t=0, calculating second-order central differences u ss , u φφ  with respect to s, φ, calculating a gradient u t , using u t =2(u ss +u φφ )−λ(u−u 0 ), where λ is a constant parameter, and updating the random event data using u(s, φ, t 2 )=u(s, φ, t 1 )+Δt u t , where Δt=t 2 −t 1 , t 1 =0 in a first iteration, and Δt is greater than 0. The method repeats the steps of calculating the second-order central differences, calculating the gradient, and updating the random event data until a change in u(s, φ, t) from a previous iteration is less than a predetermined threshold value.

BACKGROUND

1. Field

The present disclosure generally relates to random estimation inPositron Emission Tomography (PET) reconstruction.

2. Background

The most accurate and commonly used method for random estimation inPositron Emission Tomography (PET) reconstruction is the delayedcoincidence window estimation, in which the correlation between thepaired photons generated from a single annihilation can be totallyremoved, thus leaving only random events. However, due to the limitedcoincidence window and short acquisition time, statistically, themeasured delay data is only one realization of the true random datadistribution. Therefore, if not processed properly, the variance willincrease in the prompt data after subtracting the delayed coincidencewindow data.

There are many image smoothing techniques developed in a uniformlysampled data space, such as Fourier analysis, to remove high frequencynoise. The directly acquired delay raw sinogram, however, is not in auniformly sampled space. If the delay raw sinogram could be convertedinto a uniformly sampled interpolated-sinogram, then various standardsmoothing techniques could be used. The data could then beback-interpolated to the irregular-sampled raw sinogram space afterrandom sinogram smoothing in the uniformly interpolated-sinogram space.

Noise represents the high frequency components in the spectrum of thenoisy random sinogram. Random sinogram smoothing can be applied byfiltering out the high frequencies in the Fourier domain. But thesmoothing will also cause large changes in the signals that are presentin the random sinogram distribution. Furthermore, the Fourier method isa global representation of the signal and will fail for images with anirregular mask, as in the case of random sinogram smoothing in which theinterpolated random sinogram is restricted inside a mask containing allthe measurable line-of-responses. The mask represents the measurableregion in the image. Thus, any region outside the mask will not bemeasurable and thus have pixel intensity of zero.

BRIEF DESCRIPTION OF THE DRAWINGS

The disclosure will be better understood from reading the descriptionwhich follows and from examining the accompanying figures. These figuresare provided solely as non-limiting examples of the embodiments,wherein:

FIG. 1 illustrates interpolated random sinograms before and after randomvariance reduction;

FIG. 2 illustrates profiles of φ=64 in FIG. 1, before and after randomvariance reduction;

FIG. 3 is a flow diagram of a method of an embodiment of the presentdisclosure;

FIG. 4 is a flow diagram of a method of an embodiment of the presentdisclosure; and

FIG. 5 illustrates a computer system upon which embodiments of thepresent disclosure may be implemented.

DETAILED DESCRIPTION

In one embodiment, there is provided a method and apparatus forsmoothing random event data obtained from a Positron Emission Tomography(PET) scanner. In one embodiment, the method includes obtaining initialrandom event data u(s, φ, t=0)=u₀(s, φ), corresponding to t=0,calculating second-order central differences u_(ss), u_(φφ) with respectto s, φ, calculating a gradient u_(t), usingu_(t)=2(u_(ss)+u_(φφ))−λ(u−u₀), where λ is a constant parameter;updating the random event data using u(s, φ, t₂)=u(s, φ, t₁)+Δt u_(t),where Δt=t₂−t₁, t₁=0 in a first iteration, and Δt is greater than 0, andrepeating the steps of calculating the second-order central differences,calculating the gradient, and updating the random event data until achange in u(s, φ, t) from a previous iteration is less than apredetermined threshold value.

Further, in an embodiment of the method, λ has a value between 0.05 and2.0. On the other hand, in an embodiment of the method, λ has a valuebetween 0 and 0.01. Moreover, in the method, Δt has a value between 0.05and 0.1, and u(s, φ, t) converges when

$\sum\limits_{alls}{\sum\limits_{{all}\; \phi}\left( {{u\left( {s,\phi,t_{k}} \right)} - {u\left( {s,\phi,t_{k + 1}} \right)}} \right)^{2}}$

is less than the predetermined threshold value.

In one embodiment, the method includes obtaining initial random eventdata u(s, φ, z, θ, t=0)=u₀(s, φ, z, θ), corresponding to t=0,calculating second-order central differences u_(ss), u_(φφ, u) _(zz),u_(θθ), with respect to s, φ, z, θ, calculating a gradient u_(t), usingu_(t)=2(u_(ss)+u_(φφ)+u_(zz)+u_(θθ))−λ(u−u₀), where λ is a constantparameter, updating the random event data using u(s, φ, z, θ, t₂)=u(s,φ, z, θ, t₁)+Δt u_(t), where Δt=t₂−t₁, t₁=0 in a first iteration, and Δtis greater than 0, and repeating the steps of calculating thesecond-order central differences, calculating the gradient, and updatingthe random event data until a change in u(s, φ, z, θ, t) from a previousiteration is less than a predetermined threshold value.

Further, in an embodiment of the method, λ has a value between 0.05 and2.0. On the other hand, in an embodiment of the method, λ has a valuebetween 0 and 0.01. Moreover, in the method, Δt has a value between 0.05and 0.1, and u(s, φ, z, θ, t) converges when

$\sum\limits_{alls}{\sum\limits_{{all}\; \phi}{\sum\limits_{{all}\; z}{\sum\limits_{{all}\; \theta}\left( {{u\left( {s,\phi,z,\theta,t_{k}} \right)} - {u\left( {s,\phi,z,\theta,t_{k + 1}} \right)}} \right)^{2}}}}$

is less than the predetermined threshold value. In addition, the methodfurther includes applying a mask so that the obtained initial randomevent data is restricted to a region containing all measurableline-of-responses of the PET scanner.

The present disclosure describes, in one embodiment, a random sinogramsmoothing method for iteratively minimizing the L-2 norm of the sinogramimage gradients. The sinogram image gradient is calculated as thesummation of local gradients at each pixel in the sinogram space.

In two-dimensional (2D) smoothing, let u₀(s, φ) denote the pixelintensity value of a noisy image with s, φεΩ and let u(s, φ) denote thetrue noiseless image. The imaging model is then written as follows:

u₀(s, φ)=u(s, φ)+n(s, φ) in which u(s, φ) is additive noise (e.g.,Gaussian) and E[u₀(s, φ)]=u(s, φ).

The objective is to estimate u(s, φ) from u₀(s, φ) under the criterionof minimal L2 norm of image gradients. For clarity of the derivationbelow, u, u₀ are used to represent u(s, φ) and u₀(s, φ). Theoptimization problem is

$\begin{matrix}{{\min_{u}J} = {\min_{u}\left\{ {{\int_{\Omega}{\left( {u_{s}^{2} + u_{\phi}^{s}} \right)\ {s}{\phi}}} + {\lambda \frac{1}{2}{{u - u_{0}}}_{2}^{2}}} \right\}}} & (1)\end{matrix}$

with ∥u−u₀∥₂ ²=(u−u₀)², u_(s)=∂u/∂s, u_(φ)=∂u/∂φ.

Thus, the objective function to be minimized is the following:

J=∫_(Ω) L(x,y,u,u _(s) ,u _(φ))dsdφ=∫ _(Ω) {u _(s) ² +u _(φ) ²+λ½(u−u₀)² }dsdφ.  (2)

By using the Euler-Lagrange equation, we have

$\begin{matrix}{{u_{t} = {{2\left( {u_{ss} + u_{\phi\phi}} \right)} - {\lambda \left( {u - u_{0}} \right)}}}{where}{{u_{ss} = \frac{\partial^{2}u}{\partial s^{2}}},{u_{\phi\phi} = {\frac{\partial^{2}u}{\partial\phi^{2}}.}}}} & (3)\end{matrix}$

The variable t is introduced as an artificial time to parameterize thedescent direction, and u_(t) represents the gradient along the tdirection. u₀ represents the evolving image at t=0, which is theoriginal noisy image u(s, φ, 0)=u₀(s, φ). In the above-equations, λ is aconstraint parameter used to add a different amount of constraint whenthe smoothed image is much different from the original noisy image. Forthe case of strong noise (for example, where the noise overwhelms thesignal (or image)), λ=0 or λ is set to a very small value, for example,approximately less than 0.01, and greater than 0. In the case of weaknoise, λ is set between, approximately 0.5 and 2.

The image is then updated with the following equation:

u(s,φ,t ₂)=u(s,φ,t ₁)+(t ₂ −t ₁)u _(t).  (4)

In equation (4), Δt=t₂−t₁ represents the stability limit, which is setto a small value between, approximately, 0.05 and 0.1. As t increases, ade-noised image is obtained.

FIG. 3 is a flowchart of a method for random variance reduction of eventdata obtained from a PET scanner.

In step 100, initial random event data u(s, φ, t=0)=u₀(s, φ), isobtained at t=0.

In step 110, second-order central differences u_(ss), u_(φφ) arecalculated with respect to s, φ. Specifically, the second-order centraldifferences are calculated as follows:

u _(xx)=(u(x+1,y)−u(x,y))−(u(x,y)−u(x−1,y))

u _(yy)=(u(x,y+1)−u(x,y))−(u(x,y)−u(x,y−1)).

Next, in step 120, a gradient u_(t) is calculated usingu_(t)=2(u_(ss)+u_(φφ))−λ(u−u₀), where λ is a constant parameter.

In step 130, the random event data is updated using u(s, φ, t₂)=u(s, φ,t₁)+Δt u_(t), where Δt=t₂−t₁, where t₁=0 in a first iteration, and whereΔt is greater than 0.

In step 140, it is decided whether u(s, φ, t) converges, i.e., whetherthe change in u(s, φ) using a suitable norm, is less than a thresholdvalue. If u(s, φ, t) converges, then the process ends. However, if u(s,φ, t) does not converge, steps 110-130 are repeated. The convergence ofu(s, φ, t) can be determined by checking whether

$\sum\limits_{alls}{\sum\limits_{{all}\; \phi}\left( {{u\left( {s,\phi,t_{k}} \right)} - {u\left( {s,\phi,t_{k + 1}} \right)}} \right)^{2}}$

is less than a predetermined threshold.

FIGS. 1 and 2 show example results implementing the above method of FIG.3 on test data. FIG. 1 illustrates interpolated random sinograms beforeand after random variance reduction (from one (s, φ) slice). FIG. 2illustrates profiles of (φ=64 in FIG. 1, before and after randomvariance reduction.

In applying the above-discussed method to four-dimensional (4D)smoothing, the optimization problem for 4D data u(s, φ, z, θ) is asfollows:

$\begin{matrix}{{{\min_{u}J} = {\min_{u}\left\{ {{\int_{\Omega}{\left( {u_{s}^{2} + u_{\phi}^{s} + u_{z}^{2} + u_{\theta}^{2}} \right)\ {s}{\phi}{z}{\theta}}} + {\lambda \frac{1}{2}{{u - u_{0}}}_{2}^{2}}} \right\}}},} & (5)\end{matrix}$

where Ω represents the irregularly shaped region defined in the mask.The region defined in the mask does not have to be a rectangular cuboid.

The image is updated with the equation u(s, φ, z, θ, t₂)=u(s, φ, z, θ,t₁)+(t₂−t₁) u_(t), where u_(t)=2(u_(ss)+u_(φφ)+u_(zz)+u_(θθ))−λ(u−u₀).

FIG. 4 is a flowchart of a method for random variance reduction of eventdata obtained from a PET scanner.

In step 200, initial random event data u(s, φ, z, θ, t=0)=u₀(s, φ, z,θ), is obtained at t=0.

In step 210, second-order central differences u_(ss), u_(φφ), u_(zz),u_(θθ), are calculated with respect to s, φ, z, θ. Specifically, thesecond-order central differences are calculated as follows:

u _(ss)=(u(s+1,φ,z,θ)−u(s,φ,z,θ))−(u(s,φ,z,θ)−u(s−1,φ,z,θ)),

u _(φφ)=(u(s,φ+1,z,θ)−u(s,φ,z,θ))−(u(s,φ,z,θ)−u(s,φ−1,z,θ)),

u _(zz)=(u(s,φ,z+1,θ)−(s,φ,z,θ))−(u(s,φ,z,θ)−u(s,φz,θ−1))

u _(θθ)=(u(s,φ,z,θ+1)−u(s,φ,z,θ))−(u(s,φ,z,θ)−u(s,φ,z,θ−1)).

Next, in step 220, a gradient u_(t) is calculated usingu_(t)=2(u_(ss)+u_(φφ)+u_(zz)+u_(θθ))−λ(u−u₀), where λ is a constantparameter.

In step 230, the random event data is updated using u(s, φ, z, θ,t₂)=u(s, φ, z, θ, t₁)+Δt u_(t), where Δt=t₂−t₁, where t₁=0 in a firstiteration, and where Δt is greater than 0.

In step 240, it is decided whether u(s, φ, z, θ, t) converges, i.e.,whether the change in u(s, φ, z, θ) using a suitable norm, is less thana threshold value. If u(s, φ, z, θ, t) converges, then the process ends.However, if u(s, φ, z, θ, t) does not converge, steps 210-230 arerepeated. The convergence of u(s, φ, z, θ, t) can be determined bychecking whether

$\sum\limits_{alls}{\sum\limits_{{all}\; \phi}{\sum\limits_{{all}\; z}{\sum\limits_{{all}\; \theta}\left( {{u\left( {s,\phi,z,\theta,t_{k}} \right)} - {u\left( {s,\phi,z,\theta,t_{k + 1}} \right)}} \right)^{2}}}}$

is less than a predetermined threshold. The method may additionallyinclude a step of applying a mask so that the obtained initial randomevent data is restricted to a region containing all measurableline-of-responses of the PET scanner.

The aforementioned second-order partial derivaties u_(ss), u_(φφ),u_(zz), u_(θθ) are calculated on the voxels confined by the mask. If thevoxels are a boundary of the mask, the partial derivative is set to zeroin the direction where it crosses the boundary.

An advantage of the discussed embodiments is that the methods areflexible to apply any mask.

The aforementioned methods described above can be implemented using acomputer system or programmable logic. FIG. 5 illustrates a computersystem 1201 upon which embodiments of the present disclosure may beimplemented. The computer system 1201 may include, for example, aprocessing circuit, module, or unit, which perform the above-describedprocess.

The computer system 1201 includes a disk controller 1206 coupled to thebus 1202 to control one or more storage devices for storing informationand instructions, such as a magnetic hard disk 1207, and a removablemedia drive 1208 (e.g., floppy disk drive, read-only compact disc drive,read/write compact disc drive, compact disc jukebox, tape drive, andremovable magneto-optical drive). The storage devices may be added tothe computer system 1201 using an appropriate device interface (e.g.,small computer system interface (SCSI), integrated device electronics(IDE), enhanced-IDE (E-IDE), direct memory access (DMA), or ultra-DMA).

The computer system 1201 may also include special purpose logic devices(e.g., application specific integrated circuits (ASICs)) or configurablelogic devices (e.g., simple programmable logic devices (SPLDs), complexprogrammable logic devices (CPLDs), and field programmable gate arrays(FPGAs)).

The computer system 1201 may also include a display controller 1209coupled to the bus 1202 to control a display 1210, such as display 42 ora liquid crystal display (LCD), for displaying information to a computeruser. The computer system includes input devices, such as a keyboard1211 and a pointing device 1212, for interacting with a computer userand providing information to the processor 1203. The pointing device1212, for example, may be a mouse, a trackball, a finger for a touchscreen sensor, or a pointing stick for communicating directioninformation and command selections to the processor 1203 and forcontrolling cursor movement on the display 1210.

The computer system 1201 performs a portion or all of the processingsteps of the present disclosure in response to the processor 1203executing one or more sequences of one or more instructions contained ina memory, such as the main memory 1204. Such instructions may be readinto the main memory 1204 from another computer readable medium, such asa hard disk 1207 or a removable media drive 1208. One or more processorsin a multi-processing arrangement may also be employed to execute thesequences of instructions contained in main memory 1204. In alternativeembodiments, hard-wired circuitry may be used in place of or incombination with software instructions. Thus, embodiments are notlimited to any specific combination of hardware circuitry and software.

As stated above, the computer system 1201 includes at least one computerreadable medium or memory for holding instructions programmed accordingto the teachings of the present disclosure and for containing datastructures, tables, records, or other data described herein. Examples ofcomputer readable media are compact discs, hard disks, floppy disks,tape, magneto-optical disks, PROMs (EPROM, EEPROM, flash EPROM), DRAM,SRAM, SDRAM, or any other magnetic medium, compact discs (e.g., CD-ROM),or any other optical medium, punch cards, paper tape, or other physicalmedium with patterns of holes.

Stored on any one or on a combination of computer readable media, thepresent disclosure includes software for controlling the computer system1201, for driving a device or devices for implementing the invention,and for enabling the computer system 1201 to interact with a human user.Such software may include, but is not limited to, device drivers,operating systems, and applications software. Such computer readablemedia further includes the computer program product of the presentdisclosure for performing all or a portion (if processing isdistributed) of the processing performed in implementing the invention.

The computer code devices of the present embodiments may be anyinterpretable or executable code mechanism, including but not limited toscripts, interpretable programs, dynamic link libraries (DLLs), Javaclasses, and complete executable programs. Moreover, parts of theprocessing of the present embodiments may be distributed for betterperformance, reliability, and/or cost.

The term “computer readable medium” as used herein refers to anynon-transitory medium that participates in providing instructions to theprocessor 1203 for execution. A computer readable medium may take manyforms, including but not limited to, non-volatile media or volatilemedia. Non-volatile media includes, for example, optical, magneticdisks, and magneto-optical disks, such as the hard disk 1207 or theremovable media drive 1208. Volatile media includes dynamic memory, suchas the main memory 1204. Transmission media, on the contrary, includescoaxial cables, copper wire and fiber optics, including the wires thatmake up the bus 1202. Transmission media also may also take the form ofacoustic or light waves, such as those generated during radio wave andinfrared data communications.

Various forms of computer readable media may be involved in carrying outone or more sequences of one or more instructions to processor 1203 forexecution. For example, the instructions may initially be carried on amagnetic disk of a remote computer. The remote computer can load theinstructions for implementing all or a portion of the present disclosureremotely into a dynamic memory and send the instructions over atelephone line using a modem. A modem local to the computer system 1201may receive the data on the telephone line and place the data on the bus1202. The bus 1202 carries the data to the main memory 1204, from whichthe processor 1203 retrieves and executes the instructions. Theinstructions received by the main memory 1204 may optionally be storedon storage device 1207 or 1208 either before or after execution byprocessor 1203.

The computer system 1201 also includes a communication interface 1213coupled to the bus 1202. The communication interface 1213 provides atwo-way data communication coupling to a network link 1214 that isconnected to, for example, a local area network (LAN) 1215, or toanother communications network 1216 such as the Internet. For example,the communication interface 1213 may be a network interface card toattach to any packet switched LAN. As another example, the communicationinterface 1213 may be an integrated services digital network (ISDN)card. Wireless links may also be implemented. In any suchimplementation, the communication interface 1213 sends and receiveselectrical, electromagnetic or optical signals that carry digital datastreams representing various types of information.

The network link 1214 typically provides data communication through oneor more networks to other data devices. For example, the network link1214 may provide a connection to another computer through a localnetwork 1215 (e.g., a LAN) or through equipment operated by a serviceprovider, which provides communication services through a communicationsnetwork 1216. The local network 1214 and the communications network 1216use, for example, electrical, electromagnetic, or optical signals thatcarry digital data streams, and the associated physical layer (e.g., CAT5 cable, coaxial cable, optical fiber, etc.). The signals through thevarious networks and the signals on the network link 1214 and throughthe communication interface 1213, which carry the digital data to andfrom the computer system 1201 may be implemented in baseband signals, orcarrier wave based signals. The baseband signals convey the digital dataas unmodulated electrical pulses that are descriptive of a stream ofdigital data bits, where the term “bits” is to be construed broadly tomean symbol, where each symbol conveys at least one or more informationbits. The digital data may also be used to modulate a carrier wave, suchas with amplitude, phase and/or frequency shift keyed signals that arepropagated over a conductive media, or transmitted as electromagneticwaves through a propagation medium. Thus, the digital data may be sentas unmodulated baseband data through a “wired” communication channeland/or sent within a predetermined frequency band, different thanbaseband, by modulating a carrier wave. The computer system 1201 cantransmit and receive data, including program code, through thenetwork(s) 1215 and 1216, the network link 1214 and the communicationinterface 1213. Moreover, the network link 1214 may provide a connectionthrough a LAN 1215 to a mobile device 1217 such as a personal digitalassistant (PDA) laptop computer, or cellular telephone.

While certain embodiments have been described, these embodiments havebeen presented by way of example only, and are not intended to limit thescope of the inventions. Indeed the novel methods and systems describedherein may be embodied in a variety of other forms; furthermore, variousomissions, substitutions, and changes in the form of the methods andsystems described herein may be made without departing from the spiritof the inventions. The accompanying claims and their equivalents areintended to cover such forms or modifications as would fall within thescope and spirit of the inventions.

1. A method for smoothing random event data obtained from a PositronEmission Tomography (PET) scanner, the method comprising: obtaininginitial random event data u(s, φ, t=0)=u₀(s, φ), corresponding to t=0;calculating second-order central differences u_(ss), u_(φφ) with respectto s, φ; calculating a gradient u_(t), usingu_(t)=2(u_(ss)+u_(φφ))−λ(u−u₀), where λ is a constant parameter;updating the random event data using u(s, φ, t₂)=u(s, φ, t₁)+Δt u_(t),where Δt=t₂−t₁, t₁=0 in a first iteration, and Δt is greater than 0; andrepeating the steps of calculating the second-order central differences,calculating the gradient, and updating the random event data until achange in u(s, φ, t) from a previous iteration is less than apredetermined threshold value.
 2. The method according to claim 1,wherein λ has a value between 0.05 and 2.0.
 3. The method according toclaim 1, wherein λ has a value between 0 and 0.01.
 4. The methodaccording to claim 1, wherein Δt has a value between 0.05 and 0.1. 5.The method according to claim 1, wherein u(s, φ, t) converges when$\sum\limits_{alls}{\sum\limits_{{all}\; \phi}\left( {{u\left( {s,\phi,t_{k}} \right)} - {u\left( {s,\phi,t_{k + 1}} \right)}} \right)^{2}}$less than the predetermined threshold value.
 6. An apparatus forsmoothing random event data obtained from a Positron Emission Tomography(PET) scanner, the apparatus comprising: a processing circuit configuredto obtain initial random event data u(s, φ, t=0)=u₀(s, φ), correspondingto t=0, calculate second-order central differences u_(ss), u_(φφ), withrespect to s, φ; calculate a gradient u_(t), usingu_(t)=2(u_(ss)+u_(φφ))−λ(u−u₀), where λ is a constant parameter, updatethe random event data using u(s, φ, t₂)=u(s, φ, t₁)+Δt u_(t), whereΔt=t₂−t₁, t₁=0 in a first iteration, and Δt is greater than 0, andrepeat the calculating the second-order central differences, calculatingthe gradient, and updating the random event data until a change in u(s,φ, t) from a previous iteration is less than a predetermined thresholdvalue.
 7. The apparatus according to claim 6, wherein λ has a valuebetween 0.05 and 2.0.
 8. The apparatus according to claim 6, wherein λhas a value between 0 and 0.01.
 9. The apparatus according to claim 6,wherein Δt has a value between 0.05 and 0.1.
 10. The apparatus accordingto claim 6, wherein u(s, φ, t) converges when$\sum\limits_{alls}{\sum\limits_{{all}\; \phi}\left( {{u\left( {s,\phi,t_{k}} \right)} - {u\left( {s,\phi,t_{k + 1}} \right)}} \right)^{2}}$is less than the predetermined threshold value.
 11. A method forsmoothing random event data obtained from a Positron Emission Tomography(PET) scanner, the method comprising: obtaining initial random eventdata u(s, φ, z, θ, t=0)=u₀(s, φ, z, θ), corresponding to t=0;calculating second-order central differences u_(ss), u_(φφ), u_(zz),u_(θθ), with respect to s, φ, z, θ; calculating a gradient u_(t), usingu_(t)=2(u_(ss)+u_(φφ)+u_(zz)+u_(θθ))−λ(u−u₀), where λ is a constantparameter; updating the random event data using u(s, φ, z, θ, t₂)=u(s,φ, z, θ, t₁)+Δt u_(t), where Δt=t₂−t₁, t₁=0 in a first iteration, and Δtis greater than 0; and repeating the steps of calculating thesecond-order central differences, calculating the gradient, and updatingthe random event data until a change in u(s, φ, z, θ, t) from a previousiteration is less than a predetermined threshold value.
 12. The methodaccording to claim 11, wherein λ has a value between 0.05 and 2.0. 13.The method according to claim 11, wherein λ has a value between 0 and0.01.
 14. The method according to claim 11, wherein Δt has a valuebetween 0.05 and 0.1.
 15. The method according to claim 11, wherein u(s,φ, z, θ, t) converges when$\sum\limits_{alls}{\sum\limits_{{all}\; \phi}{\sum\limits_{{all}\; z}{\sum\limits_{{all}\; \theta}\left( {{u\left( {s,\phi,z,\theta,t_{k}} \right)} - {u\left( {s,\phi,z,\theta,t_{k + 1}} \right)}} \right)^{2}}}}$is less than the predetermined threshold value.
 16. The method accordingto claim 11, further comprising: applying a mask so that the obtainedinitial random event data is restricted to a region containing allmeasurable line-of-responses of the PET scanner.
 17. An apparatus forsmoothing random event data obtained from a Positron Emission Tomography(PET) scanner, the method comprising: a processing circuit configured toobtain initial random event data u(s, φ, z, θ, t=0)=u₀(s, φ, z, θ),corresponding to t=0; calculate second-order central differences u_(ss),u_(φφ), u_(zz), u_(θθ), with respect to s, φ, z, θ; calculate a gradientu_(t), using u_(t)=2(u_(ss)+u_(φφ)+u_(zz)+u_(θθ))−λ(u−u₀), where λ is aconstant parameter; update the random event data using u(s, φ, z, θ,t₂)=u(s, φ, z, θ, t₁)+Δt u_(t), where Δt=t₂−t₁, t₁=0 in a firstiteration, and Δt is greater than 0; and repeat the calculating thesecond-order central differences, calculating the gradient, and updatingthe random event data until a change in u(s, φ, z, θ, t) from a previousiteration is less than a predetermined threshold value.
 18. Theapparatus according to claim 17, wherein λ has a value between 0.05 and2.0.
 19. The apparatus according to claim 17, wherein λ has a valuebetween 0 and 0.01.
 20. The apparatus according to claim 17, wherein Δthas a value between 0.05 and 0.1.
 21. The apparatus according to claim17, wherein u(s, φ, z, θ, t) converges when$\sum\limits_{alls}{\sum\limits_{{all}\; \phi}{\sum\limits_{{all}\; z}{\sum\limits_{{all}\; \theta}\left( {{u\left( {s,\phi,z,\theta,t_{k}} \right)} - {u\left( {s,\phi,z,\theta,t_{k + 1}} \right)}} \right)^{2}}}}$is less than the predetermined threshold value.
 22. The apparatusaccording to claim 17, wherein the processing circuit is furtherconfigured to apply a mask so that the obtained initial random eventdata is restricted to a region containing all measurableline-of-responses of the PET scanner.